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The  work  documented  in  this  report  is  part  of  a  study 
of  the  travel-time  and  amplitude  anomalies  observed  at  the 
Large  Aperture  Seismic  Array  (LASA)  in  Montana.  Addenda 
will  be  made  to  this  report  as  more  data  becomes  available 
for  analysis. 

The  work  was  performed  by  the  Applied  Research  Section, 
Earth  Sciences,  a  Teledyne  Company,  316  Montgomery  Street, 
Alexandria,  Virginia,  under  Lincoln  Laboratory  Contract 
No.  BB-246 . 

This  report  was  written  by  P.  W.  Broome,  F.  A. 
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project  director  was  P.  W.  Broome. 


ABSTRACT 


Average  amplitude  anomalies  are  presented  as  functions 
of  source  region.  Data  from  approximately  220  events  were 
used  in  this  analysis. 


accepted  for  the  Mr  Force 

i’rahklin  C.  Hudson 

jhief ,  Lincoln  Laboratory  Office 


INTRODUCTION 


1  . 

Observed  short  period  teleseismic  P-phase  amplitudes 
vary  greatly,  even  between  seismometers  that  are  relatively 
closely  spaced.  These  variations  exceed  anything  that  can  be 
easily  attributed  to  instrumentation,  radiation  pattern  ef¬ 
fects,  or  obvious  geologic  factors.  These  variations  have 
been  termed  "amplitude  anomalies".  By  observing  a  large 
number  of  events,  certain  consistencies  are  found,  in  par¬ 
ticular,  some  stations  exhibit  a  preference  for  events  from 
a  given  source  region. 

The  object  of  this  report  is  to  define  the  LASA  ampli¬ 
tude  anomalies  as  functions  of  distance  and  azimuth  to  the 
source  region.  Statistical  qualifications  on  the  observed 
anomalies  and  on  the  data  from  approximately  220  events 
(Figure  1)  which  were  used  to  derive  these  anomalies  are  also 
given . 

2 .  PROCEDURE 

2 . 1  Data  Collection 

Figure  2  shows  a  typical  record  of  an  event  as 
seen  on  the  LASA  16mm.  Develocorder  film.  Two  such  films  re¬ 
cord  the  earth  motion  at  the  center  seismometer  of  each  of 
the  21  LASA  subarrays.  The  maximum  P-wave  amplitude  and  period 
within  the  first  3  or  4  seconds  of  the  signal  were  measured 
and  converted  to  ground  motion* .  Since  all  sensors  do  not  re- 

-  1  - 


*The  maximum  P-wave  amplitude  is  defined  to  be  that  phase  which 
has  the  largest  amplitude,  on  the  average,  for  all  recording 
seismometers . 
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Figure  2.  Typical  Record  Showing  Amplitude  Measurements 


cord  with  the  same  signal-to-noise  characteristics,  each  trace 
was  given  a  quality  grade  by  the  analyst  when  read. 

Quality  Grade  1  =  a  perfect  trace. 

Quality  Grade  2  =  a  trace  on  which  the  signal  is  easily 

read  and  has  a  similar  shape  (signature) 
to  the  other  traces. 

Quality  Grade  3  =  a  trace  on  which  the  signal  is  seen 

but  it  is  not  distinct  and/or  has  a 
markedly  dissimilar  signature  to  the 
other  traces. 

Quality  Grade  4  =  a  trace  on  which  the  signal  can  not  be 

seen . 

In  this  report  all  but  Grade  4  data  were  used. 

These  data  were  then  keypunched  on  IBM  cards  as  shown  in 
Figure  3  for  use  in  further  processing. 

2 . 2  Data  Reduction 

In  order  to  compare  data  from  different  events, 

it  was  first  necessary  to  remove  the  effect  of  the  sources  on 

the  received  amplitudes.  In  doing  this,  we  have  assumed  that 

the  aperture  of  LASA  is  small  compared  to  any  lobe  structure 

in  the  radiation  pattern,  and  thus  that  all  seismometers  were 

sampling  the  same  radiation  field.  The  normalization  method 

which  we  used  was  to  normalize  the  amplitude  at  each  station 

for  each  event  to  the  geometric  mean  of  the  recorded  amplitudes 

for  that  event.  Thus  if  Y, .  and  L  are  the  normalized  and 

ID  1 D 

observed  amplitudes  for  event  i  and  station  j  respectively  and 
N  is  the  total  number  of  observed  amplitudes  for  event  i,  then 
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The  derivation  and  rationalization  of  this  normalization  method 
is  given  in  Appendix  A. 

A  computer  program  was  used  to  normalize  amplitudes 
for  each  event* .  A  sample  output  from  this  program  is  shown 
in  Figure  4.  (Additional  data  computed  in  the  program  is  also 
shown ) 

2.3  Data  Display 

The  events  were  grouped  into  cells  by  distance  and 
azimuth.  With  the  assumption  that  the  amplitude  anomaly 
within  each  cell  is  constant,  a  computer  program  was  used  to 
compute  the  geometric  mean  of  the  normalized  amplitudes  for 

each  station**.  The  program  also  computed  a  measure  of  the 

♦ 

deviation  about  the  estimated  mean  as  well  as  the  number  of 
events  recorded  for  each  station  and  cell***. 

Figure  5  shows  a  sample  of  this  information  as  the 
computer  output  listing  with  appropriate  descriptions  added. 

The  data  for  all  cells  and  subarrays  are  presented  in  Tables  1 
through  21. 


-  3  - 


* 

** 


*** 


A  description  of  program  ’ANOMALY"  is  given  in  Appendix  B 
Appendix  A 

A  description  of  program  "DISPLAY"  is  given  in  Appendix  C 
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Figure  4.  Sample  Output  of  Program  “DISPLAY" 
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Figure  5.  Sample  Output  from  Program  "ANOMALY 
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302  EVENTS  PROCESSED 


3. 


ERROR  DETECTION 


Using  the  assumption  that  the  true  amplitude  anomalies 
are  constant  within  geographic  cells,  it  is  possible  to  ex¬ 
amine  the  scatter  of  the  normalized  amplitudes  in  each  cell 
for  the  purpose  of  detecting  errors  due  to  blunders  in  data 
handling  (and  other  possible  sources  of  error) .  To  detect 
and  eliminate  the  scatter  due  to  these  human  errors,  the  error 
model  study  described  in  Appendix  D  was  undertaken.  The  re¬ 
sults  of  this  study  indicate  that  the  signal  model  is  essential¬ 
ly  of  the  form  L. .  =  S.g.e. .,  where  L. ,  is  the  observed 

ID  i  D  ID  ID 

amplitude  for  event  i  and  station  j,  is  the  radiated  signal 

of  event  i,  g.  is  the  gain  factor  of  station  j,  and  e.  .  is  an 
D  ID 

error  term  which  has  a  log-normal  distribution. 


The  procedure  used  to  detect  errors  was  as  follows: 

2 

The  estimates  of  the  variances,  S  ,  computed  for  each  geo¬ 
graphic  cell  and  station  were  plotted  on  a  histogram  as  a 
function  of  the  number  of  events,  N,  in  that  cell.  By 

T 

pooling  these  variances  a  very  good  estimate,  S  ,  of  the 

2 

variance  of  the  error  population,  a  ,  was  found.  Confidence 

2  2 
intervals  for  the  S  were  then  computed  and  all  S  which 

fell  outside  of  the  upper  95%  confidence  limit  were  questioned. 

The  amplitude  measurements  involved  in  computing  these 

outlier-estimates  were  re-examined  for  blunders.  (Readings 

which  did  not  clearly  include  blunders  were  not  changed,  so 

this  re-examination  of  the  data  should  not  bias  the  final 

results.)  The  percentage  of  reading  blunders  so  detected  was 


4 


less  than  1%  of  the  total  readings  taken.  With  these  errors 
eliminated,  the  entire  process  was  repeated  with  the  results 
shown  by  Figure  6.  On  this  second  pass  the  total  number  of 
outliers  above  the  upper  confidence  limit  is  2.1%  of  the 
total  number  of  estimates.  This  is  in  good  agreement  with 
the  idealized  distribution  (i.  e.,  2.1%  versus  2.5%). 

4 .  CONCLUSIONS 

It  is  clear  from  the  results  of  this  study  that  certain 
stations  at  the  Large  Aperture  Seismic  Array  show  a  preference 
"on  the  average"  for  events  from  particular  regions.  Due  to 
the  high  variability  of  the  data,  however,  one  might  hypothe¬ 
size  that  the  averaging  yields  only  the  low  frequency  components 
of  an  anomaly  versus  distance-azimuth  function  which  contains, 
in  addition,  a  considerable  amount  of  high  frequency  components. 
The  amount  of  available  data  is  not  adequate  to  provide  an 
answer  to  this  hypothesis.  It  is  interesting  to  note,  however, 
that  the  average  anomalies  for  the  eight  Fiji  events,  used  in 
Appendix  A  to  establish  the  distribution  of  the  errors,  have 
the  same  distribution  function  as  the  errors;  see  Figures  A3 
and  A4 . 

The  average  anomalies  are  large  in  many  cases,  with 
factors  of  3  to  1  between  stations  being  common.  As  such  they 
must  play  a  part  in  studies  involving  the  estimation  of  source 
radiation  patterns.  The  variability  of  the  data  might  also 
put  limits  on  the  degree  to  which  one  can  hope  to  estimate  these 
radiation  patterns  using  any  seismic  net:  a  trade-off  between 
resolution  and  accuracy  likely  exists. 
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Figure  6.  Distribution  of  Standard  Deviation  Estimates  Computed 

from  N  Events 


The  signal  model  which  was  derived  for  error  detection 
served  well  in  that  function.  The  model  is  also  interesting 
in  its  own  right,  however.  The  physical  significance  of  the 
multiplicative  error  factor  was  not  to  be  part  of  this  study, 
and  for  this  reason  no  effort  was  made  to  explain  it.  It  is 
possible  that  the  effect  could  have  some  bearing  upon  the 
placement  of  new  arrays. 

The  data,  as  presented,  is  adequate  to  describe  the 
average  amplitude  anomalies  in  the  more  seismically  active 
regions . 

As  more  events  are  analyzed,  it  will  be  possible  to  re¬ 
duce  the  variability  of  the  average  amplitude  anomaly  estimates, 
especially  in  those  cells  for  which  only  a  small  number  of 
events  were  available.  Since  the  average  anomaly  varies  from 
cell  to  cell,  subdivisions  within  a  cell  may  exhibit  dis¬ 
similar  average  anomalies.  Thus,  as  more  data  is  accumulated, 
cells  might  necessarily  be  re-defined  in  order  to  obtain  opti¬ 
mum  cell  distribution  consistent  with  event  distribution. 
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APPENDIX  A 


A . 1  INTRODUCTION 

The  signal  levels  recorded  at  the  stations,  aside  from 
statistical  variations,  will  be 


L.  .  =  S  g. 
13  1  3 


( A-l ) 


where  is  the  signal  radiated  by  the  ith  event  and  g_.  is 

a  gain  associated  with  the  jth  station.  The  set  of  ratios 
of  these  levels  for  a  single  event  are  termed  "amplitude 
anomalies" .  We  define  an  individual  anomaly  for  stations  t 
and  m  as 


A ,  =  L  „  .  /L  . 

<?,mi  ti  mi 


(A-2) 


note  that 


A  .  =  (A,  . ) 

m£i  •cmi 


-1 


(A-3) 


and  denote  the  set  of  these  as 


A 

1 


.  ;A  \  21 
l  ^miJ-?.  ,m=l 


(A-4) 


All  elements  in  A^  are  equally  important  but  only  20 

of  the  441  are  independent.  Any  20  which  contain  each  station 
index  at  least  once  are  sufficient  to  construct  the  entire 
set . 


A-l 


In  analyzing  the  recorded  signals,  it  was  found  that 

was  highly  variable.  To  reduce  this  variability,  averages 

were  taken  over  events  from  the  same  region.  This  pooling  of 
events  is  possible  under  the  assumption  that  amplitude  anoma¬ 
lies  are  a  slowly  varying  function  of  distance  and  azimuth* . 
Thus,  if  the  region  is  sufficiently  small,  this  is  equivalent 
to  assuming  the  anomalies  are  constant  for  each  region. 

Several  approaches  were  taken  to  this  pooling  problem. 
Each  of  these  was  equivalent  in  the  sense  that  they  all  pro¬ 
duce  the  same  final  results.  A  summary  of  these  will  be 
presented  in  the  remainder  of  this  section. 

A. 2.1  Method  1 

From  Equation  (A- 3)  we  note  that  the  apparent 
size  of  an  anomaly  depends  upon  the  reference  station.  To 
eliminate  this  effect,  a  logarithmic  transformation  was  made. 

Since,  on  a  logarithmic  scale,  the  reference 
station  which  is  chosen  serves  only  to  translate  all  anoma¬ 
lies  equally,  we  can  dispense  with  a  reference  station  momen¬ 
tarily  and  concentrate  instead  on  a  set  of  translating  factors. 

First,  consider  the  problem  of  fitting  the  ob¬ 
served  amplitudes  for  a  single  event  to  a  selected  level  S. 


A-2 


*This  assumption  is  open  to  question.  It  is  clear  from  the 
results  of  this  study,  however,  that  certain  stations  do  show 
a  preference  "on  the  average"  for  events  from  particular 
regions.  From  the  data,  though,  one  might  hypothesize  that 
the  averages  yield  only  the  slowly  varying  components  of  a 
distance-azimuth  function  which  contains  a  considerable  amount 
of  high  frequency  components.  The  amount  of  available  data 
is  not  adequate  to  answer  this  question. 


We  seek  a  multiplier  or  normalizing  factor,  ,  to  minimize 
the  total  deviation  from  S, 


e .  =  Y  (log  S-log  G.L. .)2  (A-5) 

1  L>  1  lj 

j 

stations 


and  find  that 


G. 

1 


J 

=  S/  (  17 

j=l 


1/J 


V 


(A-6) 


where  J  is  the  total  number  of  stations. 

Second,  consider  minimizing  the  total  vari¬ 
ation  between  events.  We  define  this  total  variation  as 

I  Uog  G.L..  -  log  GkLk.)  (A-7) 

j  i,k 

station  events 


where  the  G  are  multiplicative  factors  as  before. 

This  minimum  will  be  found  among  the  solutions 
to  the  equations 


5  V 

BG 

n 


0, 


n  =  1,2,... 


(A-8) 


or 


I 


2 

G 

n 


( log  G 
r  n 


iog  Gk  Lk.) 


0 


(A-9) 


A- 3 


which  reduces  to 


I  log  GnJ  =  £  log 

j 


Gk  Lkj 


n: 


(A-10) 


where  I  and  J  are  the  total  number  of  events  and  stations 
respectively. 


By  setting 

£  log  Gk  I*  =  IJ  log  S, 
j  ,k 

a  constant,  we  obtain 


(A-ll) 


1/J 


G  =  S/  (  7T  L .  .  ) 

1  j-i  13 


(A-12) 


as  before. 


The  constant  S  can  be  chosen  arbitrarily.  The 


weights  G^  are  independent  and  the  normalizing  process  is 

equivalent  to  fitting,  in  a  least-square  sense,  each  set  of 

measurements,  L. .,  to  a  reference  level  Log  S.  This  permits 

13 

addition  or  omission  of  events  from  a  collection  without  dis¬ 
turbing  the  normalization  factors  on  the  remaining  data.  In 
this  way  we  obtain  a  set  of  normalizing  factors  which  can  be 
used  on  an  event  by  event  basis. 

The  normalized  amplitudes  can  be  plotted  as 
ordinates  on  semi-log  graph  paper  with  stations  as  abscissae 
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on  the  linear  axis.  The  resultant  plots  will  be  adjusted  in 
a  least-squares  sense  to  minimize  the  total  variation. 
Outliers  will  tend  to  be  made  more  evident  on  such  plots. 

To  use  these  data  for  computing  anomalies  we 
first  define  the  average  normalized  amplitude  at  each  station 
as 


I 

log  g.  =  <  log  G .  L . .  >  =  —  Y  log  G.  L. .  (A-13) 

j  ii]  I  Li  ii] 

i=l 

events 


where  the  G^  are  previously  determined  normalizing  factors. 
Then  the  average  anomaly  is  given  as 


<  A  .  .  >  =  log  f  <  log  G .  L .  >  -  <  log  G .  L . ,  > 

mh  L  1  lm  i  il  J 

1/1  (A-14) 


I 

=  TT 

i=l 


L 

lm 
Li  l 


which  is  the  geometric  mean  of  the  individual  anomaly  esti¬ 
mates.  This  reasonable  result,  along  with  the  desirable 
features  of  the  logarithmic  transformation  are  reason  enough 
for  using  this  method. 

A. 2. 2  Method  2 

Up  to  this  point  in  the  analysis  it  was  not  neces¬ 
sary  to  incorporate  an  error  model.  To  lend  further  just  if i- 
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cation  to  the  above  method,  however,  it  was  felt  desirable  to 
do  this.  We  therefore  assumed  an  error  model  of  the  form 


L.=S  g  .  e  .  .  +  e  .  . 
ID  1  d  ID  ID 


(A-15) 


where  e.  .  is  due  to  the  additive  factors  of  noise  and 
ID 

reading  errors  and  e  accounts  for  variations  in  system  gain 

due  to  equipment  changes,  effects  of  wave-shape,  etc. 

The  steps  taken  will  be  to  assume  and  then 
justify  two  hypotheses  which  were  put  forth  after  looking  at 
a  large  amount  of  data.  These  assumptions  are  that  the  addi¬ 
tive  noise  e. .  is  negligible  for  large  events  and  that  e. 

ID  ID 

has  a  log  normal  distribution. 

By  neglecting  the  additive  error,  one  may 
linearize  the  model  to 


log  L^_.  =  log  +  log  g_.  +  log  e^_. 


(A-16) 


and  by  minimizing  the  sum  of  squares 


L.  . 


=  I  dog  =  l  (1o9  ) 


1 >  D 


1/  D 


!g3 


( A-17 ) 


estimate  the  values  of  and  g_.  (albeit  not  independently). 

These  critical  values  of  gain  and  signal  level 
are  in  the  solutions  to  the  equations 


SE 

as 


=  0  and 


n 


Be 


m 


=0  n=l,2,...,I, 

m  =  1 , 2  ,  .  .  .  ,  J 


( A-18 ) 


A-6 


or 


y  Log  (L  /S  g  )  =  0 
L-t  nj  n  j 


(A-19.1) 


and 


Y  Log  (L  /S  g  )  =  0 
L->  im  1  m 


(A-19.2) 


which  can  be  rewritten  as 


1/J 


S  =  (  rr  L  .  /  g.) 
n  j=l  D 


(A-20.1) 


and 


m 


1/1 


g  =  (  n  L  /  s. ) 


i=l 


im  l 


(A-20.2) 


Then  by  adding  one  constraint,  namely  that 


TT  s 


j=l 


we  have 


(A-21) 


1/1 


gm  = 

m 


(S  tt  im) 
i=l 


IJ 


T/KT 


( A-22 ) 


(  tt  L.  .  ) 

ij-i  13 
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which  is  equivalent  to  Equation  (A-12)  when  1=1  and 

g  =  G  L , 
ro  1  1m 

The  anomalies  are,  as  before,  ratios  of  the 
estimates  of  the  gain  factors, 


<  Am-C,j  >  ^m  /  ^ l 


1/1 


n 

i=l 


(L.  /  L  ) 

lm  il 


(A-2  3) 


( A-24 , 
A-14) 


A . 3  Distribution  of  Errors 

To  establish  the  distribution  function  for  e. .  , 

3-D 

eight  events  from  the  Fiji  Islands  were  normalized  using 
(A-6)  with  S  =  1  as  shown  in  Figure  A1  .  The  average  gain 
of  each  station  was  computed  using  (A-22),  and  the  deviations 
from  these  averages  were  plotted  as  shown  in  Figure  A2  . 

This  data  was  then  broken  into  class  intervals  and  plotted  on 
a  logarithmic  probability  graph,  Figure  A3  . 

The  fit  to  a  log-normal  distribution  function  is  suf¬ 
ficiently  good  to  justify  the  assumptions  that  the  additive 
noise  is  negligible  for  our  purposes  provided  we  use  large 
events.  Also,  since  the  error  term  in  the  linearized  model 
is  distributed  normally,  classical  least  squares  methods  with 
the  usual  statistical  qualifications  on  all  estimates  can  be 
applied  directly. 

Another  way  of  presenting  this  data  is  shown  in 
Figure  A4  .  There  the  observed  amplitudes  for  these  events 
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are  plotted  against  the  model  values,  namely  observed  log  L^_. 


versus  log  +  log  g^  where  and  g^  are  estimated  using 
Equations  (A-20)  and  (A-21)  - 
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NORMALIZED  AMPLITUDES 


STATION 


Figure  A-l .  Eight  Fiji  Island  Events  -  Normalized  and  Showing  the 

Means  for  Each  Station 
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Figure  A-2  .  Eight  Fiji  Island  Events  Showing  Deviations  about  the 

Mean  Amplitudes 
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Figure  A-4 .  Observed  versus  Predicted  Signal  Levels  -  3  Fiji  Island 
Events.  Using  a  Product  Error  Model  with  Zero  Additive  Noise  Assumed. 
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APPENDIX  B 


SEISMIC  DATA  LABORATORY 
ALEXANDRIA,  VIRGINIA 

DIGITAL  COMPUTING  SECT I ON 


A.  IDENTIFICATION 

Title;  ANOMALY 

COOP  Identification;  Z108 

Category:  General 

Programer :  Nicholas  Fletcher 

Date:  16  June  1966 


B .  PURPOSE 


To  compute  travel-time  and  amplitude  anomalies 
C .  USAGE 

1.  Operational  Procedure:  This  is  a  FORTRAN  63  program. 

All  input  is  from  punched  cards  and  the  output  is 
from  the  printer,  and  an  off-line  card  image  tape. 

2.  Parameters :  Input  parameters  consist  of  the  latitude 

and  longitude  of  the  location  of  the  subarray 
centers . 


3. 

Space  Required:  Not  Applicable 

4. 

Temporary  Storage  Requirements: 

None 

5. 

Alarms:  None 

6  . 

Error  Codes:  None 

7  . 

Error  Stops:  None 

8. 

Input  and  Output  Tape  Mounting: 

There  is  no  input 

tape.  The  output  tape  is  on 

logical  Unit  16. 

9. 

Input  and  Output  Formats: 

B-l 


The  input  card  formats  are  as  follows: 

Subarray  locations  deck,  must  be  the  first  deck 
of  cards  read  by  the  program.  The  first  card  of 
this  deck  must  contain  the  number  of  subarray 
location  cards.  This  number  is  punched  in 
Columns  1  through  5  right  justified.  The  first 
card  is  followed  by  the  subarray  location  cards, 
with  one  location  per  card.  The  program  is 
limited  to  a  maximum  of  51  subarrays.  1  (One) 
punched  in  Column  10,  prints  out  J.B.  Table  and 
station  locations. 

Subarray  location  card  format  is  as  follows: 

FORMAT  (A8,  2  (F5.0,  F3.0,  F5.1,  Al) 

COLS. 

1-8  Subarray  name,  left  justified 

9-13  Subarray  latitude  degrees 

14-16  Subarray  latitude  minutes 

17-21  Subarray  latitude  seconds  (with  decimal  point 

punched) 

22  N  or  S  for  north  or  south  latitude 

23-27  Subarray  longitude  degrees 

28-30  Subarray  longitude  minutes 

31-35  Subarray  longitude  seconds  (with  decimal  point 

punched) 

36  E  or  W  for  east  or  west  longitude 

The  location  cards  are  followed  by  the  first 
event  card  - 


Event  card  format  is  as  follows: 

Columns  1  through  14  contain  the  origin  time  of 
the  event 

COLS. 

1-2  Month 
Day 


3-4 
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COLS. 


5-6  Year 

7-8  Hour 

9-10  Minute 

11-14  Second  to  the  nearest  tenth  of  a  second 

15-18  Latitude  of  the  event  in  degrees  to  the  nearest 

tenth 

19  N  or  S  for  north  or  south 

20-24  Longitude  of  the  event  in  degrees  to  the  nearest 
tenth 

25  E  or  W  for  eat  or  west  longitude 

26-28  Depth  of  event  to  the  nearest  kilometer 

29-31  Magnitude  to  the  nearest  tenth 

32-36  This  is  an  output  field.  The  program  will  com¬ 
pute  the  distance  (to  the  nearest  kilometer) 
between  the  event  and  the  standard  subarray 
(AO  for  LASA) 

37-39  This  is  an  output  field.  The  program  will  com¬ 
pute  the  station  to  epicenter  azimuth  (to  the 
nearest  degree  relative  to  AO) . 

40-78  Not  used  by  the  program,  but  reproduced  on 
output . 

79-80  Must  be  blank,  to  identify  this  card  as  an  event 
card  and  not  a  phase  card . 

v1* 

Each  event  card  should  be  followed  by  phase 
cards.  A  phase  card  contains  the  observed  phase 
data  for  a  given  subarray.  Each  phase  card 
must  have  an  index  number  punched  in  Columns 
79-80.  This  index  is  established  by  the  order 
of  the  subarray  location  cards  (See  Sec-C9)  . 
There  does  not  have  to  be  a  phase  card  for 
every  subarray  location  card  nor  a  specific 
numerical  order. 


Phase  card  format: 

FORMAT  (Columns  1  through  9  contain  the  phase  arrival 
time  at  this  subarray) 
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COLS. 


1-2  Hour 

3-4  Minute 

5-9  Second  to  the  nearest  hundredth 

10  Not  Used 

11-14  Amplitude  (zero  to  peak)  to  the  nearest 

millimicron 

15-16  Not  Used 

17-20  Amplitude  (peak  to  peak)  to  the  nearest 
millimicron 

21  Not  Used 

30  Signal  Quality 

31-35  Time  anomaly  to  the  nearest  hundredth  of  a 

second — this  is  an  output  item  computed  by  the 
program. 

36-40  Amplitude  anomaly  zero  to  peak  relative  to 

geometric  mean  computed  by  the  program  (output 
field) 

41-42  Subarray  name  code  (output  item) 

43  Not  Used 

44-45  Reference  subarray  (for  each  event)  (only  in 

last  phase  card) 

46-50  Peak  to  peak  amplitude  anomaly  relative  to 
geometric  mean  (output  field) 

51-55  Zero  to  peak  amplitude  anomaly  normalized  to 
average  amplitude  (0-P)  (output  field) 

56-60  Peak  to  peak  amplitude  anomaly  normalized  to 
average  amplitude  (P-P)  (output  field) 

61-78  Reserved  but  not  used  by  program 

79-80  Index  number  of  subarray  recording  this  phase 

must  correspond  to  subarray  location  card 

Additional  events  with  their  associated  phase 
cards  follow  in  like  fashion  for  as  many  events 
as  desired.  The  subarray  location  cards  are 
not  repeated  in  the  deck  setup. 
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The  output  tape  is  a  card  image  tape,  with 
the  same  information  and  format  as  the  input 
tape  with  the  addition  of  the  computed  values 
for  distance  and  azimuth  for  event  cards,  and 
travel-time  and  amplitude  anomalies  for  the 
phase  cards. 

10.  Selective  Jump  and  Stop  Settings;  None 

11.  Timing:  Approximately  3.5  seconds  per  event 

12.  Accuracy:  No  better  than  the  input  data 

13.  Cautions  to  User:  Every  phase  card  must  have  an  index 

number.  If  Columns  79-80  are  blank,  the  program 
assumes  it  is  an  event  card 

14.  Equipment  Configuration:  Standard  F-63  COOP  monitor 

system  with  the  card  reader  on  Unit  7  and  the 
printer  on  Unit  6  and  the  output  tape  on  Unit  16 

15.  References :  Letter  from  MIT  dated  April  8,  1966 


D.  METHOD 


Definitions 


Observed  Travel-Time  =  Observed  Phase  Arrival  Time  Minus 

Event  Origin  Time 

Predicted  Travel-Time  =  Travel-Time  for  the  J-B  Table. 

Found  by  table  look-up  on  distance  and  depth  with 
four-point  interpolation  for  distance  and  linear 
interpolation  for  depth. 

0^  =  Observed  travel-time  at  subarray  i 
=  Predicted  travel-time  to  subarray  i 

0^  =  Observed  travel-time  at  standard  subarray 
=  Predicted  travel-time  to  standard  subarray 

Transit  Time  Residual  (Travel-Time  Anomaly)  =0.-0  +H  -H. 

l  r  r  i 

(Relative)  Amplitude  Anomaly  =  Amplitude  (Zero  to  Peak)  at 
subarray  i  divided  by  geometric  mean  amplitude  (zero  to 
peak) 


B-5 


(Average)  Amplitude  Anomaly  =  Amplitude  (Zero  to  Peak)  at 
subarray  i  divided  by  arithmetic  average  amplitude 
(zero  to  peak)  for  all  functioning  subarrays. 

Relative  and  average  amplitude  anomalies  are  also  computed 
in  like  manner  for  the  peak  to  peak  measurements. 
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SEISMIC  DATA  LABORATORY 
ALEXANDRIA,  VIRGINIA 

DIGITAL  COMPUTING  SECTION 


A.  IDENTIFICATION 

Title;  Display  and  Analysis  of  Anomalies 
COOP  Identification:  Z95  DISPLAY 
Category:  General 

Programer:  J.  W.  Monroe 

Date :  19  July  1966 

B .  PURPOSE 

To  take  the  output  from  the  anomaly  program  and  com¬ 
pute  the  average,  standard  deviation,  and  number  of 
occurrences  of  a  selected  anomaly  for  selected  subarrays 
or  sensors,  distances,  and  azimuths. 

C .  USAGE 


1.  Operational  Procedure?  This  FORTRAN-63  program  is 

with  a  binary  or  symbolic  deck  and  a  set  of 
control  cards  for  each  use. 

2 .  Parameters :  None 

3.  Space  Required:  21,222  buffers 

4.  Temporary  Storage  Requirement:  None 

5.  Print-Out s :  The  program  first  prints  out  the  anomaly 

used,  followed  by  the  distance  and  azimuth  criteria 
and  the  subarrays  or  sensors  used.  If  any  seismo¬ 
grams  are  to  be  deleted  from  processing,  they  are 
listed  next.  Following  this  will  be  a  page  for 
each  subarray  or  sensor  with  the  average,  standard 
deviation  and  number  of  occurrences  of  the  selected 
anomaly  for  each  distance  and  azimuth  chosen.  If 
plotting  is  requested,  the  last  page  will  contain 
scale  factors  and  other  information  pertaining  to 
the  plots. 
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Error  Returns:  None 


7.  Error  Stops:  None 

8 .  Input  and  Output  Tape  Mountings: 

Logical  Unit  #4  -  Input  BCD  tape  from  ANOMALY 
Logical  Unit  #3  -  Output  plot  tape  if  desired 

9 .  Input  and  Output  Formats: 


INPUT: 

(only 

one  case 

per  run) 

Card  # 

Cols 

Format 

Description 

1 

1 

11 

Anomaly  to  be  used  (see  last  page) 

5 

11 

Plot  switch  1-plot,  blank-no  plot 

2 

1-10 

110 

Number  of  Dels  (Max  19) 

11-20 

110 

Number  of  Azimuths  (Max  13) 

21-30 

no 

Number  of  Subarrays  or  Sensors 

(Max  21) 

next  N 

1-10 

F10.2 

Distance  (km) 

11-20 

F10.2 

+_  Distance  (km) 

Where  N 

=  NO. 

punched 

in 

Cols  1-10  on  card  #2  (Max  19) 

next  M 

1-10 

F10.2 

Azimuth  (Deg) 

11-20 

F10.2 

+_  Azimuth  (Deg) 

Where  M 

=  NO. 

punched 

in 

Cols  11-20  on  card  #2  (Max  13) 

next  L 

1-10 

F10.0 

Subarray  or  Sensor  Number 

Where  L 

=  NO. 

punched 

in 

Cols  21-30  on  card  #2  (Max  21) 

next  D 

1-4 

14 

Seismogram  Nos.  to  be  deleted 

(Max  100) 
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10. 

11. 

12  . 
13. 


14. 

15. 


OUTPUT;  The  printed  output  will  be  that  described  in 
PRINT -OUTS .  The  plot  tape  on  logical  unit  #3  will 
contain  a  plot  for  each  subarray  or  sensor  of  the 
average  for  that  subarray  minus  the  average  of 
averages  for  all  subarrays.  The  plotted  value  will 
be  computed  for  each  distance  at  each  azimuth  (in 
that  order)  until  N  points  have  been  plotted  where 
N  equals  the  number  of  dels  times  the  number  of 
azimuths.  The  plots  will  always  represent  a  new 
del  every  tenth  of  an  inch. 

Selective  Jump  and  Stop  Settings:  None 

Timing:  2  mili-sec  per  case  where  each  case  has  a 

separate  distance,  azimuth,  subarray,  and  anomaly. 

Accuracy:  The  same  as  input  data 

Cautions  to  Users:  Both  logical  units  3  and  4  must 
be  mentioned  in  the  COOP  card.  If  no  plot  is 
desired,  use  the  BY3  convention.  The  plot  routine 
is  not  internally  controlled  so  the  user  must  re¬ 
quest  the  same  number  of  files  to  be  plotted  as 
there  are  subarrays  or  sensors  with  a  maximum  of 
ten  across  the  plot.  Observe  that  there  are  maxi¬ 
mum  limits  on  the  number  of  distance  and  azimuth 
intervals  and  number  of  subarrays  and  seismograms 
to  be  deleted. 

Equipment  Configuration:  Standard  COOP  for  FORTRAN-63 

References : 

Program  write-up  for  ANOMALY  Fletcher  16  June  1966 


10. 

11. 

12  . 
13. 


14. 

15. 


D .  METHOD 


X. 


Standard  Deviation  =  a 


i ,  j  ,k 


N  .  -  1 


Number  of  Occurrences  =  N 


i ,  j  ,k 
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where,  i  =  azimuth  +_  deg's 
j  =*  distance  +  km's 
k  =  subarray  or  sensor 

N  =  number  events  pertinent  to  the  specified  coordinates 
and  specified  subarray 

X.  .  =  anomaly  used 

i »  ]  i  k 
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SAMPLE  INPUT  DATA  FOR  DISPLAY 


S'  EOT 


O 

l 

<_n 


Event  nos.  punched  on 
cards  (1  per  card)  shall 
not  be  processed  Max  (100) 


Subarray  No . 


£  80.0  20.0 
s'  150.0  10.0 

' 200.0  5070 

315.0  15.0  Format  (2F10.2) 


AZI  Intervals 

315  +  15  Deg 


etc . 


10000.0  200.0 

8000.0  200.0 


DEL  Intervals 
Format  (2F10.2)  8000  +_  200  km 


Si  4 


1 

Execute , 


5  •* -  Format  (3110) 


-  Format  (ll,3X,Il) 

50,  56,  1. 


3  =  Number  of  Del  Intervals  -  Max (19) 

4  =  Number  of  Azi  Intervals  -  Max (13) 

5  =  Number  of  Subarrays  or  Sensors  -  Max (21) 

Anomaly  To  Be  Used 

1  -  Time  Plot  Switch  (col  no.  5) 

2  -  Amplitude  (Z  -  P)  Relative  blank  -  no  plot 

3  -  "  (P  -  P)  "  1  -  plot  desired 

4  -  "  (Z  -  P)  Average 

5  -  "  (P  -  P) 


TO: 

FROM: 

DATE : 
SUBJECT  s 
REFERENCE: 


Users  of  Program  DISPLAY 
J.  W„  Monroe 
2  September  1966 

Additional  Features  of  Program  DISPLAY 
Z95  DISPLAY  Program  Writeup,  19  July  1966 


The  user  may  now  run  DISPLAY  using  more  than  one  input  tape  de¬ 
scribed  in  the  program  writeup  referenced  above.  The  first  input  tape 
must  be  on  Logical  Unit  4  but  the  second,  third,  etc.,  tapes  may  be 
assigned  to  any  logical  unit  from  1  thru  49  inclusive,  except  Logical 
Unit  3  (the  plot  tape).  When  using  more  than  one  input  tape,  Logical 
Units  3,  4,  and  all  other  logical  units  must  be  mentioned  in  the  COOP 
card.  The  control  card  setup  for  the  case  of  more  than  one  input  tape 
is  exactly  the  same  as  that  shown  in  the  SAMPLE  INPUT  DATA  FOR  DISPLAY 
of  the  writeup,  except  that  following  the  last  card  of  seismograms  to 
be  deleted  or  following  the  last  subarray  number  if  there  are  no  de¬ 
letion  cards,  a  blank  card  must  be  inserted  followed  by  the  logical 
tape  assignment  for  the  second  input  tape.  If  more  than  two  input 
tapes  are  used,  additional  cards  are  used  to  specify  the  logical  tape 
unit  for  input  tapes  3,  4,  etc.,  but  the  blank  card  is  not  used  agciin. 
Format  for  these  tape  assignments  is  15.  The  following  is  a  sample 
card  setup  starting  with  the  deletion  cards: 


etc 


£ 


£ 


£ 


fourth  input  tape 


logical  tape  unit  for  the 
third  input  tape 


logical  tape  unit  for  the 
second  input  tape 


r 


c 


r 


|  BLANK  CARD  «= — 

seismograms  to  be  deleted 

'  1350 

4210 


3490 
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As  soon  as  an  end  of  file  is  encountered  on  the  input  taps  being 
processed,  a  tape  assignment  card  is  read  in  and  the  procram  pauses 
with  the  next  input  logical  unit  number  displayed  in  the  "A"  register. 
When  the  Start  Key  is  hit,  the  program  begins  processing  the  new  input 
tape.  Two  input  tape  units  might  be  used,  such  as  4  anu  5,  where  the 
first  tape  is  on  Unit  4  and  the  second  tape  is  mounted  sn  5  while  4  is 
being  processed.  After  Unit  5  begins  processing,  anotler  tape  could  be 
mounted  on  Unit  4. 


DEBUGGING  AIDS  FOR  DISPLAY 


Jump  Key  2  On  displays  the  event  number  being  processed  in  octal 

in  the  "A"  register.  Program  halts  until  Start  Key 
is  hit.  For  exact  stops,  leave  JK-2  or  and  step 
through  events  with  Start  Key.  Program  halts  on  header 
card  of  event  being  displayed  in  "A"  register. 


Jump  Key 


Jump  Key 


Up  causes  printout  of  accumulation  thc.t  far  into  the 
program.  Turn  JK-1  off  as  soon  as  printout  begins 
and  JK-3  off  as  soon  as  printout  has  finished.  Use 
this  feature  sparingly  since  some  loss  of  accuracy 
occurs  when  sigma  is  set  equal  to  zero  rather  than  a 
negative  number.  All  other  times,  the  accumulation 
is  reset  to  its  original  value  without  loss  of  accu¬ 
racy. 


Be  sure  both  COOP  card  and  EXECUTE  card  have  sufficient  time  limits  for 
the  case  being  run.  For  a  case  of  324  events  with  21  stations,  17 
distances,  and  12  azimuths,  allow  at  least  60  minutes. 


CHANGE  IN  FORMAT 


The  format  of  the  control  cards  containing  the  numbers  of  seismo¬ 
grams  to  be  deleted  has  be  changed  from  14  to  15. 
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APPENDIX  D 


I 


4 


*• 


The  product  error  model  described  in  Appendix  A  can  be 
used  to  detect  amplitude  data  which  is  inconsistent  with  the 
population  from  which  it  is  drawn.  Serveral  approaches  to 
automating  this  error  detection  were  suggested  and  the  follow¬ 
ing  is  a  description  of  the  method  which  was  finally  selected 
for  production  use. 


The  collection  of  recorded  data  was  grouped  into  a  set 
of  cells.  The  grouping  geographically  was  by  distance  and 
azimuth  and  a  cell  is  defined  by  a  particular  distance  and  azi¬ 
muth  range  and  a  station  designation.  The  assumptions  made 
were  that  the  geographic  regions  chosen  are  sufficiently  small 
so  that  the  true  anomalies  are  constant  within  each  region 
and  that  the  variability  is  the  same  for  all  cells. 

Using  the  results  of  Appendix  A,  namely  that 


log  e  . 
ID 


N  (0,a  ),  one  may  make  estimates,  denoted  as  S  , 


of  a  for  each  cell.  These  can  be  pooled  to  obtain  a  very 

2 

good  estimate  of  a  ,  and  with  this  estimate  and  known 

distribution,  outlying  estimates  and  thus  probable  blunders 
can  be  detected. 


The  estimates  of  a  ,  computed  for  each  cell  (given 
station  j  and  region)  are 


I 


1 


D  — 1 


where  I  is  the  number  of  events  recorded  in  the  particular 
cell , 


X. 

1 


log 


J =1 


(  TT  L.  .) 

ID 

j=l 


1/J 


X. 

l 


2  2 
and  S  is  distributed  as  Xj  ^ 


By  pooling  these  variance  estimates,  a  very  good 
overall  estimate  is  obtained  due  to  the  large  number  of  data. 
The  estimates  in  each  geographical  cell  were  quantized, 
indexed  and  entered  on  a  histogram,  (Figure  6  of  the  main 
report)  as  a  function  of  the  number  of  events  in  that  cell. 

A  pooled  estimate  of  the  population  variance  is  then 
given  by 


GO  GO 

-ll  l  * 

1=2  J=0 


1  »  j 


00 

I 


1=2  J=0 


where  f  .  is  the  number  of  estimates  in  the  class  interval 
I/D 

indexed  by  j  and  with  I-events  in  the  cells.  The  value  of 
2 

S .  is  equal  to  the  center  value  of  the  variance  in  the 
D 

class  interval,  j. 


D-2 


2 

Using  this  estimate  for  a  we  may  compute  100(l-a)% 


ii* 


. 


2  2  2 
confidence  intervals  for  S„  and  S„  .  Thus,  we  set  a  =  S 

N  N 

and  define  confidence  intervals 


Ct)  * 


(a/2,f ) 


<  SN2  <  Q  )  K2 


f  (l-n/2,£) 


or 


(ft/2, f) 


Q  <  SN  < 


*  (l-a/2,f) 
f 


H 


where  f  =  number  of  degrees  of  freedom,  (I  -  1),  and 


X  and  x  *  are  found  in  tables 

<«/2,f>  ( 1-/2 , f ) 
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